//----------------------------------*-C++-*----------------------------------//
/*!
 * \file   MomentsIndex.hh
 * \author Jeremy Roberts
 * \date   Aug 30, 2011
 * \brief  MomentsIndex class definition.
 * \note   Copyright (C) 2011 Jeremy Roberts.  Inspired by LM_Indexer of
 *         Denovo.
 */
//---------------------------------------------------------------------------//
// $Rev:: 132                                           $:Rev of last commit
// $Author:: j.alyn.roberts@gmail.com                   $:Author of last commit
// $Date:: 2011-09-07 19:49:16 +0000 (Wed, 07 Sep 2011) $:Date of last commit
//---------------------------------------------------------------------------//


#ifndef MOMENTSINDEX_HH_
#define MOMENTSINDEX_HH_

#include <vector>

#include "tran/Definitions.hh"
#include "tran/Traits.hh"

#include "../utilities/DBC.hh"
#include "../utilities/SP.hh"

namespace slabtran
{

//===========================================================================//
/*!
 *  \class MomentsIndex
 *  \brief Indexing for moments field access.
 *
 *  The one group (or within-group) scattering source in 3-d is defined
 *  \f[
 *      Q(\mathbf{r},\mathbf{\Omega}) =
 *      \sum^L_{l=0} \frac{2l+1}{4\pi}
 *      \sum^l_{m=-l} \Sigma_s{l}(\mathbf{r}) \phi^m_l(\mathbf{r}) \: .
 *  \f]
 *  Such sums over \f$l\f$ and \f$m\f$ are commonplace, in generating
 *  sources and also in simply accessing moments sequentially.  This class
 *  builds a vector of \f$(l,m)\f$ pairs that can be used in a single loop
 *  over all flux moments \f$\phi^m_l\f$.
 *
 *  For 3-d, the moments are ordered as
 *  \f[
 *      [(0,0)] \, , \, \, \, [(1,-1),\, (1,0)\, (1,1)] \, ,
 *      \, \, \, [(2,-2),\,(2,-1)\, \ldots \: .
 *  \f]
 *  For a Legendre order of \f$L\f$, there are \f$(L+1)^2\f$ moments.
 *
 *  For 2-d, several moments can be eliminated by symmetry.  Only those
 *  moments for which \f$l+m\f$ is even are retained (see Hebert).  Thus,
 *  the moments are ordered as
 *  \f[
 *      [(0,0)] \, , \, \, \, [(1,-1),\, (1,1)] \, ,
 *      \, \, \, [(2,-2),\,(2,0)\, \ldots \: .
 *  \f]
 *  For a Legendre order of \f$L\f$, there are \f$(L+1)(L+2)/2\f$ moments.
 *
 *  For 1-d, \f$m=0\f$, so that the moments are ordered by \f$l\f$ alone,
 *  from 0 to the given order.  Hence, there are \f$L+1\f$ moments.
 *
 *  \tparam     D   Dimension.
 */
template <class D>
class MomentsIndex
{

public:

  /// Typedefs
  //\{
  //
  typedef typename util::SP<MomentsIndex<D> >   SP_momentsindex;
  //\}

  /*!
   *  \brief Constructor.
   *
   *  \param legendre_order     Legendre order of the \em flux.
   */
  explicit MomentsIndex(int legendre_order);

  /// \name Accessors
  //\{

  /*!
   *  \brief Return vector of values of m for a given l.
   *
   *  \param  l  Legendre degree.
   *  \return    m values
   */
  const std::vector<int>& m_index(int l) const;

  /*!
   *  \brief Return legendre order.
   */
  int legendre_order() const
  {
    return d_legendre_order;
  }

  /*!
   *  \brief Return l value.
   *  \param    index   Cardinal moment index.
   *  \return           l value.
   */
  const int l(int index) const
  {
     Require(index >= 0);
     Require(index < d_number_moments);
     return d_l[index];
  }

  /*!
   *  \brief Return m value.
   *  \param    index   Cardinal moment index.
   *  \return           m value.
   */
  const int m(int index) const
  {
    return d_m[index];
  }

  // \}

private:

  /// Legendre order of the \em flux.
  const int d_legendre_order;

  /// Number of moments.
  const int d_number_moments;

  /// Vector of proper \f$m\f$ values for a given \f$l\f$.
  std::vector<std::vector<int> > d_m_index;

  /// Vector of l values.
  std::vector<int> d_l;

  /// Vector of m values.
  std::vector<int> d_m;


};

} // end namespace slabtran

//---------------------------------------------------------------------------//
// INLINE FUNCTIONS
//---------------------------------------------------------------------------//

#include "angle/MomentsIndex.i.hh"

#endif /* MOMENTSINDEX_HH_ */

//---------------------------------------------------------------------------//
//              end of MomentsIndex.hh
//---------------------------------------------------------------------------//
